SIM202 : HOME SAVER¶

Cyril NEFZAOUI & Bertrand PATUREL

Présentation des codes et des résultats

Sommaire :¶

Visualisation des données

Imputation par continuité

Modèle linéaire

GAM

GBM

XGBoost

Imputation VS prévision pure

Permutation

Complétion partielle

Initialisation¶

On récupère les deux jeux de données train (jeu de données d'entrainement) et test (jeu de données de submission).

Visualisation des données¶

On regarde tout d'abord a quoi ressemble notre jeu de données Appliances en fonction du temps.

Puis l'on s'interresse au jeu de donnée test. On souhaite comprendre quelle est la plage des relevés de cette table et la comparer aux dates du jeu de données train.

Conclusion : On suppose donc qu'il y avait un jeu de données complet initialement qui à été répartie en deux jeux de données, l'un train et l'autre test

Avant tout traitement de train et test, nous avons rajouté ces data frame d'une colonne $is in test$, dont les coefficients valent 1 ou 0 selon que l'élèment associé appartient au train ou au test. On a ensuite fusionné train et test dans une large $data frame$ que l'on appelle total. On effectue un tri de total par ordre chronologique croissant, ce qui nous permet de visualiser la répartion de certaines valeurs du test dans le train.

Traitement des données¶

On cherche en priorité à savoir ou sont les données manquantes dans les jeu de données.

On remarque que les colonnes RH 6 et Visibility qui correspondent à la visibilité et à l'humidité sur le facade extérieure de la maison ont de valeurs manquantes. Affichons ces deux variables pour comprendre pourquoi il y a des données manquantes :

On comprend donc que les valeurs manquantes sont dus à la saturation des capteurs de l'humidité et de la visibilité lorsqu'il y a des conditions extrêmes.

On choisit alors deux méthodes pour compléter les valeurs saturées.

• Première méthode : Pour _RH6 on fixe simplement à 10 les valeurs manquantes
• Deuxième méthode : Pour Visibility on prend la moyenne pondérée aux bords non saturés : $V_t^{NA}=\frac{k}{j+k}V_{t-j}+\frac{j}{j+k}V_{t+k}$ où $(V_t)$ représente la série des visibilités.

Imputation statistique des valeurs maquantes par continuité¶

Nous venons de mettre en évidence qu'il y avait deux parties, une partie comprenant l'imputation statistique et l'autre de la prédiction. Dans cette partie nous nous concentrons sur une méthode concise.

Notons $(A_t)$ la série de toutes les Appliances, $(A_t^e)$ celle du test à compléter et $(A_t^r)$ celle du train. Pour chaque valeur Appliances $(A_t^e)$ manquante dans le test, on regarde les deux valeurs de train $A_{t+k}^r$ et $A_{t-j}^r$ encadrant et on effectue une moyenne pondérée de ces deux valeurs. Puis on améliore cette aproximation en attribuant un poids différent à ces deux Appliances encadrantes. Ainsi plus la valeur manquante à compléter du test est proche temporellement de son bord gauche (resp droit) plus le poids attribuer au bord gauche (resp droit) est important.

On utilise donc la relation suivante : $A_t^e$=$\frac{k}{j+k}A_{t-j}^r+\frac{j}{k+j}A_{t+k}^r$

Le modèle linéaire¶

Choix du modèle linéaire¶

Ici on teste tout d'abord plusieurs manière d'effectuer un STEPAIC. On choisit finalement la méthode forward qui nous laisse plus aisément le choix des variables.

On cherche désormais à ajuster le nombre de degrés de libertés du modèle de la régression linéaire. On se propose tout d'abord un petit code pour voir le nombre de degrés de liberté minimisant le rmse sur un échantillon proche des valeurs à prédire.

On cherche donc à calculer la trace de la matrice $X(XX^{\top})^{-1}X^{\top}$. On obtient alors 12 degrés de libertés, l'on choisit donc de prendre 12 varaiables explicatives.

On remarque que les données ont du mal à absorber les pics. La regression linéaire à du mal à prédire les valeurs extrêmes.

Remarque : On voit clairement que notre prédiction peut être améliorée

Finalement on a un score Kaggle - RMSE : 67.9

Amélioration du modèle linéaire par classification¶

Une fois notre meilleur modèle linéaire choisit on décide de l'améliorer en classifiant les données. Nous montrons plus tard l'importance de certaines données, notamment la température de la cuisine et son humidité (qui indiquent clairement que la famille est en train de préparer le repas par exemple) où la température et l'humidité de la salle de bain. Si ces varaiables ont des valeurs fortes nous sommes sur des pics pour appliances et ainsi nous boostons les pics de notre régression linéaire à l'aide d'une méthode proche de celle des K plus proche voisins.

Notre code comporte 4 étapes :

• si les habitants dorment, le modèle prend des valeurs de consommation faibles;
• si les habitants cuisinent, on choisit de placer un pic;
• si les habitants se lavent, on choisit de placer un pic;
• si il ya du repassagen on place un pic;
• si les habitants partent en week end, le modèle prend des valeurs de consommation plus faibles.

Finalement on a un score Kaggle - RMSE : 67.8

Modèle GAM¶

Choix du modèle GAM¶

Nous avons choisi de faire un modèle GAM. Nous souhaitons savoir quelle variable explicative mettre dans le modèle et confimer nos intuitions.

On décide alors de choisir en partie NSM et lights comme variables explicatives (bien que deux variables corrélés n'implique que l'une explique l'autre) pour commencer. Nous avons ensuite complété notre choix de variables avec un StepGam. Voici en noir la courbe reélle d'une partie des valeurs de test ayant servi d'entrainement et en rouge la prédiction de ses valeurs par le gam.

Finalement on a un score Kaggle - RMSE : 70

Modèle mixte : GAM et Linéaire¶

Nous avons améliorer nos deux modèles de prédiction GAM et linéaire en faisant une combinaison de ces mmodèles. L'idée est de laisser le modèle GAM s'exprimer lorsqu'il y a des pics importants et lorsque les habitants dorment (pour combler les défaults).

On effectue donc tout d'abord une prédiction $(g_t)$ avec un modèle gam expliquant les Appliances uniquement par NSM. On souhaite avoir un modèle prédictif valant zéro lorsque les habitants dorment et valant 1 lorsqu'il y a un pic de consommation. On pose donc $G_t=\frac{g_t-min(g_t)}{max(g_t)-min(g_t)}$.

Nos Appliances prédites seront alors $(A_t)$ telque : $A_t=51*(1-G_t)+L_t*G_t$

où $L_t$ correspond à notre meilleur prédiction linéaire des Appliances manquantes du test.

Justification de ce modèle : Lorsque les habitants dorment, $G_t$ est proche de $0$ donc $A_t$ équivaut à $51$ ce qui est en moyenne la valeur des Appliances de minuit à 7h00 du matin. Quand il y a un pic $G_t$ est proche de 1 et donc $A_t$ correspon à la prédiction linéaire $L_t$.

Finalement on a un score Kaggle - RMSE : 67.9

Utilisation d'un GBM¶

On peut visualiser un indicateurs de la performance : l'erreur Out Of Bag qui doit normalement converger vers 0 quand on augmente le nombre d'itération dans l'algorithme GBM.

Definition : L'Erreur Of Bag est $\mathbb{E}((m(X)-Y)^2)$ ou $m(X)$ est la prédiction en fonction des variables explicatives $X$ et $Y$ la variable cible (ici Appliances). L'Erreur OOB représente l'écart moyen entre la prédiction et la variable cible, en régression on choisit naturellement pour formule $\frac{1}{n}\sum_{k=1}^n(\hat{Y}_k-Y_k)^2$ avec $\hat{Y}=m(X)$.

#On réduit le jeu de données, on garde celle qui nous interresse et on enleve les redondances :
library(dplyr)
train_gbm=(train %>% select(-c(date,Day_of_week,DayType,lights,Visibility,InstantF,Posan,Month,Heure,BE_load_actual_entsoe_transparency,
                               BE_load_forecast_entsoe_transparency,BE_wind_onshore_capacity,BE_wind_onshore_capacity,BE_wind_onshore_generation_actual,BE_wind_onshore_profile,WeekStatus)))



library(gbm)
object_gbm=gbm(formula =my_formula ,

               data =train_gbm, distribution = "gaussian",
               , var.monotone = NULL, n.trees = 100,
               interaction.depth = 1, n.minobsinnode = 10, shrinkage = 0.05,
               bag.fraction = 0.5, train.fraction = 0.5, cv.folds = 20,
               keep.data = TRUE, verbose = TRUE,   n.cores = 4)

On remarque bien que certaines variables ont une influence relative faible. En appliquant le principe de parcimonie, nous avons donc décidé de retirer ces variables pour ne garder que les plus significatives.

Nous remarquons une similitude : ce sont les mêmes variables qui sont retenues dans le cadre de GBM et dans le cadre du modèle linéaire obtenu par step.AIC. On peut donc entrainer un modèle linéaire obtenu par step.AIC à l'aide d'un GBM.

Utilisation de l'interface XGBoost (Xtreme Gradient Boosting)¶

Nous avons utilisé XGBoost en tant que GBM amélioré.En effet,XGBoost tout comme GBM utilisent le principe de gradient Boosting, mais il nous a semblé que XGBoost proposait une implémentation plus performante que GBM, grâce à notamment:

• parallélisation du calcul

• tous les paramètres clés de l'apprentissages sont modifiables

• utilisation de matrices pleines pour accélerer le temps de calcul.

• La modélisation a été concue de telle sorte à limiter l'overfitting

En fait, nous avons utilisé XGBoost comme une interface, qui regroupe en fait plusieurs modèles. Premièrement, on définit une structure de contrôle qui permet de choisir la méthode de validation : nous avons successivement choisi la validation croisée (cv) et l'Erreur Out of Bag (OOB).

Nous nous sommes alors demandé comment les arbres de régression développés dans la phase d'apprentissage allaient mesurer l'homogénéité dans les données de Appliances, variable quantitative.

##importation des différentes bibliothèques
library(tidyverse)
library(lubridate)
library(xts)
library(dygraphs)
library(ranger)
library(dplyr)

library(mltools)
library(data.table)
library(xgboost)
library(caret)

p2=read.csv(file="submission_73.2.csv", sep=",", dec='.')$Appliances[4655:5771]
p3=read.csv(file="submission_73.3.csv", sep=",", dec='.')$Appliances[4655:5771]
p4=read.csv(file="submission_73.9.csv", sep=",", dec='.')$Appliances[4655:5771]

p_moy=(p2+p3+p4)/3
p_moy_original=p_moy

p_moy=list(p_moy)
###Méthode de complétion partielle pour conserver les valeurs d'Appliances déjà prédites et vraisemblables
p_moy_df=setDT(p_moy, keep.rownames=FALSE, key=NULL, check.names=FALSE)
p_moy$to_use=rep(0, times = 1117)
p_moy$to_use[which(p_moy$V1<80)]=1
p_moy$V1[which(p_moy$to_use==0)]=NA
par(mfrow=c(2,2))

plot(ts(p_moy_original))
lines(ts(p_moy$V1),
      col="red")

###On cree une copie de total et on l'appelle heavy pour insister sur le fait qu'il contient toutes les variables explicatives.
total_heavy=total
total_heavy_matrix=data.matrix(total_heavy, rownames.force = NA)
total_heavy_numeric= (lapply(total_heavy, function(x) as.numeric(x)))
total_heavy_numeric_df=setDT(total_heavy_numeric, keep.rownames=FALSE, key=NULL, check.names=FALSE)

library(doParallel)
library(missForest)
registerDoParallel(cores=4)

###on creee une variable total_light qui contient uniquement les variables pertinentes.
total_light=(total_heavy_numeric_df %>% select(-c(rv1,rv2,lights,Tdewpoint,RH_6,Visibility,Instant,DayType,BE_load_actual_entsoe_transparency,
                                                                      BE_load_forecast_entsoe_transparency,BE_wind_onshore_capacity,BE_wind_onshore_capacity,BE_wind_onshore_generation_actual,BE_wind_onshore_profile,Month,Windspeed,Press_mm_hg,InstantF,T2,T3,RH_2,RH_8,Day_of_week,Posan)))
permutation=sample(1:19735)
total_light_permuted=total_light[permutation]

total_with_imputation=missForest(total_light, maxiter = 10 ,ntree = 10, variablewise = TRUE,
                                decreasing = FALSE, verbose = TRUE,

                                mtry =sqrt(ncol(total_light)),
                                replace = TRUE,
                                classwt = NULL, parallelize = 'forests')

##On recupere les donnees avec les Appliances imputees
data_complete=total_with_imputation$ximp
#######Approche Xgboost

library(xgboost) 

#####Définition de la structure de controle
xgb_trcontrol = trainControl(method = "cv", number = 1000,allowParallel = TRUE, verboseIter = TRUE, returnData = FALSE)

#####definitiion du grid
xgbGrid <- expand.grid(nrounds =1,  max_depth =50,colsample_bytree = 0.9, eta = 0.1,gamma=0,min_child_weight = 1,subsample = 0.5)

total_light$is_in_test[which(is.na(total_light$Appliances))]=2
to_train=total_light[which(total_light$is_in_test<2)]
y_train=total_light$Appliances[which(total_light$is_in_test<2)]
total_light_wt_App=(to_train %>% select(-Appliances))
X_train=(xgb.DMatrix(as.matrix(total_light_wt_App ),label=y_train))


xgb_model=xgb.train(params=list(booster="gbtree"),data=X_train, trControl = xgb_trcontrol, tuneGrid = xgbGrid,  method = "xgbTree",nrounds=100000)
rm(X_test)
X_test=total_light[which(total_light$is_in_test==2)]
X_test=(X_test %>% select (-Appliances))
X_test=xgb.DMatrix(as.matrix(X_test))
predictions=predict(xgb_model,newdata=X_test)
plot(ts(predictions))

reconstitution_test=test_total$Appliances
length(which(is.na(reconstitution_test)))
tokaggle=c(appliances_imputed,predictions)
plot(tokaggle)

total_light$Appliances[which(is.na(total_light$Appliances))]=predictions
length(which(is.na(test_total$Appliances)))
test_total$Appliances[which(is.na(test_total$Appliances))]=predictions
which(is.na(test_total$Appliances))
plot(test_total$Appliances[4655:5771])
which(is.na(total_light$Appliances))

tokaggle=test_total$Appliances[4655:5771]
plot(to_g)

Techniques de prévision hybride, imputation des valeurs manquantes vs prevision pure¶

On remarque que le test a été echantilloné de manière linéaire dans le train. On dispose donc d'un échantillonage de bonne qualité, car les modèles peuvent être entrainés sur le train, ce qui aurait été plus difficile si les valeurs manquantes de Appliances avaient été réparties sans homogénéité.

Voici le résultat de l'imputation des données manquantes à l'aide d'une méthode de prévision à l'aide de Xgboost en utilisant une méthode de validation croisée.

Par ailleurs, voici le résultat de l'imputation des données manquantes à l'aide d'une méthode de complétion par forêts aléatoires. Plus précisement, nous avons utilisé le $package missForest$. Dans cette méthode d'imputation, l'erreur Out of Bag converge vers 0, ce qui indique le validité de notre imputation.

Finalement on a un score Kaggle - RMSE : 69.1

Remarque : Nous choisissons de travailler sur les jeu de données entier : on considère que les données imputées font partie des données initiales, donc on peut les réutiliser pour entrainer nos modèles. Nous allons discuter de ce choix dans la suite de notre étude.

Tests de nos modèles à l'aide des données imputées¶

Précedemment, nous avons établi que les données imputées étaient très proches des vraies données au sens du RMSE.Dès lors, on fait l'hypothèse que ces données sont plausibles.On peut donc considérer ces données manquantes avec deux points de vue:

1) Ces données manquantes peuvent être complétées par imputation (ACP par missMDA, forêts aléatoires par missForest)

2) Ces données peuvent être prédites à l'aide d'un modèle qui a été entrainé sans prendre en compte ces données.

Sur cette figure, on a représenté 3 séries mises en jeu dans le test

• en vert, les valeurs imputées par missForest,

• en violet, la partie présente dans le train,

• en bleu,la partie de prédiction pure.

Conclusion : Pour améliorer notre score, nous décidons d'améliorer notre prédiction pure dans notre algorithme.

Nous représentons ici plusieurs prédictions, en utilisant plusieurs jeux de paramètres dans l'algorithme XGBoost

De plus, puisque tous les prédicteurs sont numériques,on se ramène en fait à une régression constante par morceaux sur des pavés de $R^{p}$,ces pavés étant déterminés par dichotomie succesives lors de la création de l'arbre de régression

Méthode de permutation-imputation¶

Nous avons constaté que, une fois le problème d'imputation résolu, le problème de prévision des valeurs dans le futur présente une difficulté supplémentaire: on doit prédire un un ensemble continu dans le temps de Appliances ,alors que dans le problème d'imputation, on pouvait compléter les valeurs manquantes grâce aux valeurs des voisins des Appliances manquantes, par exemple par la méthode des plus proches voisins.

Nous avons donc appliqué une permutation sur l'ensemble d'apprentissage afin d'homogénéiser la répartition des valeurs à prédire de facon uniforme.Ainsi, on transforme un problème de prédiction en un problème d'imputation.

En vert est représentée une prédiction réalisée à l'aide d'une régression linéaire classique.

En bleu est représentée une prédiction réalisée à l'aide de notre méthode de permutation-imputation.

Désormais, on a trois possibilités de prédire les valeurs manquantes:

1) Prévision à l'aide d'un modèle (régression, GBM, Xgboost, Random Forest)

2) Imputation, sans permutation.

3) Imputation, avec permutation.

Interprétation : Sur ce graphe, on superpose les résultats issus d'une même méthode d'imputation, avec et sans permutation. Avec permutation, les valeurs à prédire sont dispersées au sein du jeu de données et le lien temporel entre chaque relevés successif à predire est moindre. Ainsi la variance est plus élevée.

L'ACF de la série des Appliances d'entrainement montre bien une corrélation de la série à son passé.

Méthode d'amélioration des résultats : Complétion partielle¶

On dispose désormais de prévisions de bonne qualité, c'est-à-dire avec un faible RMSE. On choisit nos 3 meilleures prédictions au sens du RMSE, puis on crée la série qui correspond à la moyenne de ces prévisions. Le signal résultant est intéressant car on retrouve de faible appliances la nuit. Ainsi les valeurs d'Appliances très faibles (seuil Appliances <80) paraissent vraisemblables. Nous choisissons donc de supposer ces valeurs comme vraies. On réduit ainsi le cardinal de l'ensemble des données à prédire.En fixant ce seuil , on ne prédit par exemple que 893 valeurs au lieu de 1117. En rouge sont représentés les valeurs de Appliances qu'on introduit dans le jeu d'entrainement, et qu'on a donc plus besoin de prédire.

Conclusion¶

On peut comparer les avantages et les inconvenients de nos diverses méthodes. Dans un premier temps, nous avons mis en place une simple regression linéaire.Cette méthode a le mérite d'être facilement interprétable, et elle constitue une méthode qui conserve une continuité.En effet, de petites variations sur les données conduiront à de petites variations sur les résultats, grâce à la continuité du produit matriciel. Cependant, les régressions linéaires sont par nature inadaptées pour représenter des variations brusques (à savoir les pics).

Au contraire, les méthodes d'arbre permettent de modéliser de facon plus fidèle les discontinuités.Cependant, de faibles variation sur les données peuvent conduire à de grandes variations sur les résultats.Cependant, ces défauts ont été écartés grâce à l'utilisation du bootstrap.On a aussi établi que notre cas, l'utilisation d'un arbre de régression s'identifiait à une régréssion linéaire sur des fonctions constantes par morceaux, sur des pavés de $R^{p}$.

Remarque générale : Lors de l'utilisation de nos modèles, nous avons été confrontés au compromis biais-variance. Des modèles simples, avec peu de variables explicatives ont eu généralement un fort biais mais une variance peu élevée, et les modèles complets ont parfois réussi à obtenir un biais inférieur mais une variance élevée.

Nous sommes deuxième de la compétition.